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ABSTRACT 

We compute the structure function scaling of the integrated intensity im- 
ages of two J=l-0 13 C0 maps of Taurus and Perseus. The scaling exponents 
of the structure functions follow the velocity scaling of super-sonic turbulence, 
suggesting that turbulence plays an important role in the fragmentation of cold 
interstellar clouds. The data also allows to verify the validity of the two ba- 
sic assumptions of the hierarchical symmetry model, originally proposed for the 
derivation of the velocity structure function scaling. This shows that the same 
hierarchical symmetry holds for the projected density field of cold interstellar 
clouds. 



Subject headings: turbulence - ISM: kinematics and dynamics - radio astron- 
omy: interstellar: lines 
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1. Introduction 

Due to the complexity of the Navier-Stokes equations, mathematical work on turbulence 
is often inspired by experimental and observational measurements. Since geophysical and 
laboratory flows are predominantly incompressible, turbulence studies have been limited 
almost entirely to incompressible flows (or to infinitely compressible ones, described by the 
Burgers equation). Little attention has been paid to highly-compressible, or super-sonic 
turbulence. 

The cold interstellar medium (ISM) of galaxies such as the Milky Way is both highly 
turbulent and highly super-sonic, within a range of scales from hundreds of parsecs to 
approximately a tenth of a parsec. Large observational surveys of the cold ISM of our 
galaxy and of dust within it have become available in the last few years and new surveys of 
unprecedented sensitivity and resolution will be obtained in the near future by FIR satellites 
such as SIRTF and Herschel. 

The cold ISM provides a good laboratory for super-sonic turbulence and results from 
new observational surveys should motivate new mathematical work. Furthermore, the inter- 
pretation of the astronomical data requires a knowledge of basic properties of super-sonic 
turbulence. 

Previous works have tried to investigate the properties of ISM turbulence by estimating 
the second order structure function or the power spectrum of the velocity field in dark clouds 
(Kleiner & Dickman 1987; Hobson 1992; Miesch & Bally 1994; Brunt & Heyer 2002a,b), 
as sampled by molecular emission lines, with a number of different methods. The centroid 
velocity at each map position is used as an estimate of the local radial velocity. This centroid 
velocity results from a convolution of density, velocity and excitation temperature along the 
line of sight and is not easily related to the three dimensional velocity structure. Other 
works have instead studied the column density distribution of dark clouds, estimated from 
the integrated intensity of molecular emission lines or from dust thermal emission. Fractal 
dimensions have been estimated (Beech 1987; Bazell & Desert 1988; Scalo 1990; Dickman 
et al. 1990; Falgarone et al. 1991; Zimmermann et al. 1992; Henriksen 1991; Hetem & Lepine 
1993; Vogelaar & Wakker 1994; Elmegreen & Falgarone 1996). Multifractal (Chappell & 
Scalo 2001) and wavelet (Langer et al. 1993) analysis have also been proposed as a way of 
characterizing the projected density structure of molecular clouds. 

The multifractal analysis applied by Chappell & Scalo (2001) to dust continuum images 
is related to the scaling of the moments of the projected density. In the present work we study 
the scaling of integrated intensity differences (structure functions) of the J=l-0 transition of 
13 CO from the Taurus and the Perseus molecular cloud complexes. If we had computed the 
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moments of intensity, instead of intensity differences, then or analysis would be equivalent 
to the multifractal analysis by Chappell & Scalo (2001). 

We have already shown in previous works that super-sonic turbulence in a roughly 
isothermal gas, such as the cold ISM, generates a complex density field with density contrasts 
of several orders of magnitude and with statistical properties consistent with observational 
data from molecular clouds (Padoan et al. 1998b, 1999; Padoan & Nordlund 1999). We usu- 
ally refer to this process as turbulent fragmentation. Turbulent fragmentation of star forming 
clouds is unavoidable, since super-sonic turbulent motions are ubiquitously observed. The 
results of this work provide further evidence of the importance of turbulent fragmentation 
in star forming clouds. 

In § 2 we compute the relative scaling of the structure functions and in § 3 we discuss 
our results. The hierarchical structure model is briefly presented in § 4 and its validity for 
our data is verified with the so called (5 and 7 tests. We draw our conclusions in § 5. 



We use two observational J=l-0 13 CO spectral maps of the Perseus MC complex (Billawala 
et al. 1997) and of the Taurus MC complex (Mizuno et al. 1995). These MC complexes have 
an extension of approximately 30 pc, and rms radial velocity a v = 2.0 km/s in Perseus 
and o~ v = 1.0 km/s in Taurus. The rms sonic Mach number is therefore approximately 14 
in Perseus and 7 in Taurus, assuming a gas kinetic temperature of approximately 10 K. 
The angular resolution is approximately 0.09 pc in the Perseus map, assuming a distance of 
300 pc and 0.08 pc in the Taurus map, assuming a distance of 140 pc. 

In each map we compute an image of the integrated intensity, I(r), defined as: 



where T(v, r) is the antenna temperature at the velocity channel v and map position r. 
The integrated intensity is set to zero in each map position with peak antenna temperature 
smaller than five times the rms noise averaged over the whole map. The integrated intensity 
of J=l-0 13 CO is roughly proportional to the total gas column density since in most lines 
of sight of the two maps analyzed here this transition is not optically thick. However, the 
conversion from J=l-0 13 CO integrated intensity and gas column density also depends on 
the abundance ratio between 13 CO and H 2 , and on the precise spatial distribution of the 
excitation temperature (see Padoan et al. (1998a) for a detailed computation of the effect 
of the spatial distribution of excitation temperature on the estimation of the J=l-0 13 CO 



2. Structure Functions of Projected Density Field 
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column density). 

The structure functions of the integrated intensity image J(r) are defined as: 



^(r) = (|/(rO-/(r' + ArOn = (|A/n 



(2) 



where p is the order and the average is extended to all map positions r' and all position 
differences Ar' such that |Ar'| = r. 

It is possible to compute the moments by first deriving the probability density function 
(PDF) of AI, P(AI,r), from the observational map, and then computing: 



The PDF can be obtained by smoothing appropriately the histogram obtained from the 
observational data, as suggested by Leveque & She (1997). Here, however, we prefer to 
average the moments directly over the whole observational sample, as in (2), without first 
deriving the PDF. 

The structure functions are plotted in Figure 1, up to the order p = 20. They are well 
approximated by a power law, at least at low orders, between approximately 0.3 and 3 pc, 
for both Taurus and Perseus. We have performed a least square fit in that range of scales 
for all orders. The least square fits are shown in Figure 1 as solid lines. We call rj(p) the 
exponent of the power law fit of the p order moment 1 : 



We then compute the scaling exponents relative to the third order, ri(p)/r)(3), following the 
idea of extended self-similarity (Benzi et al. 1993; Dubrulle 1994). 

Values of r)(p)/r)(3) for both Taurus and Perseus are plotted in Figure 2, up to the 20th 
order, and are compared with those predicted for the scaling of the structure functions of the 
velocity field (((p)/((3)) by the Kolmogorov model (Kolmogorov 1941), the She & Leveque 
model (formula (7) below) and the model by Boldyrev (2002). The latter is an extension 
of the She-Leveque model to super-sonic turbulence, verified with numerical simulations of 
super-sonic turbulence by (Boldyrev et al. 2002a,b). The scaling exponents of the structure 
functions of the integrated intensity images of Taurus and Perseus follow very closely the 
velocity scaling found in Boldyrev (2002), within 5% accuracy up to the 20th order for Taurus 
and within 10% accuracy for Perseus. 



1 The Greek letters £ and t are usually reserved for the exponents of the moments of velocity differences 
and energy dissipation respectively. 
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For practical purposes it may be useful to know the intrinsic value of the scaling expo- 
nents, not just that relative to the third order one. For Taurus we find 77(3) = 1.10 and for 
Perseus 77(3) = 1.18, both close to the Kolmogorov's value for the velocity structure func- 
tion, £(3) = 1.0. The second order exponents are r](2) = 0.77 for Taurus and 77(2) = 0.83 for 
Perseus, not far from the value of 77(2) = 0.7 predicted analytically (Boldyrev et al. 2002a) 
from the velocity scaling in Boldyrev (2002). 

3. Statistical Significance of High Order Moments 

The velocity scalings of laboratory, geophysical or numerical flows are hardly verified 
up to the 10th order (Miiller & Biskamp 2000; She et al. 2001). In the present work we have 
used astronomical data to constrain moments beyond the 10th order. A way to verify the 
statistical significance of high order moments is to plot the integrand of the expression (3), 
and verify that at its peak position the PDF is defined by a significant number of samples. 
The peak position of the integrand of (3) grows with increasing order p and converges to the 
value of the cutoff of the PDF. For a finite sample size the cutoff of the PDF may depend 
on the sample size, in which case high order moments would be inaccurate. However, as 
pointed out in Leveque & She (1997), if a definite cutoff is present, for example if it is well 
defined by a significant number of samples, then the moment calculation up to any order 
should be accurate, and the fast convergence of high order moments is not an artifact of the 
finite sample size. 

We have plotted in Figure 3, left panel, the tail of the PDF of |A/| = | J(r') -I(r' + Ar')| 
from the Taurus data, for |Ar'| = 2.6 pc. Vertical dotted segments mark the peak position 
of the integrand of (3) for moments of increasing order, from p = 3 to p = 15. The moments 
converge at p — 15, and are therefore accurate up to that order. However, the PDF shows 
a very sharp cutoff. The cutoff is definitely significant since it is defined by a few hundred 
samples. The 'real' PDF describing the physical process may therefore have a similarly 
sharp cutoff, in which case the convergence around the 15th order found in the observational 
data would be significant, and orders above the 15th would be accurately estimated by the 
present data set. The moment convergence corresponds to the asymptote of the function 
i](p)/i](3) plotted in Figure 2. A similar result regarding the PDF cutoff is obtained for 
the Perseus data (right panel of Figure 3). However, in this case the peak of the integrand 
of (3) for every order is very shallow. Although formally the peak position has converged 
already at the 5th order (right panel of Figure 3), the moments converge around the 8th 
order. The convergence may again be well defined, as for the Taurus data, since the PDF 
cutoff is defined by a large sample size. 
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4. The Hierarchical Structure Model for Super-Sonic Turbulence 

The scaling of the velocity structure functions in incompressible turbulence are best 
described by the She-Leveque (1994) formula: 

C(p) 

C(3) 

where 



7P+C(l-^) (5) 



1-37 

is interpreted as the Hausdorff codimension of the support of the most singular dissipative 
structures. In incompressible turbulence the most dissipative structures are organized in 
filaments along coherent vortex tubes with Hausdorff dimension D = 1 and so C = 2. 
Furthermore, f3 3 = 2/3 (She & Leveque 1994), which yields the She-Leveque's velocity 
scaling for incompressible turbulence: 

C(v) /o\p/3 



C(3) 



(7) 



Boldyrev (2002) has proposed to apply the scaling of (5) to super-sonic turbulence, with 
the assumption that the Hausdorff dimension of the support of the most singular dissipative 
structures is D — 2 (C = 1), since dissipation of super-sonic turbulence occurs mainly in 
sheet-like shocks. Using the physical interpretation of (5) by Dubrulle (1994), 7 = 1/9 and 
C = (2/3)/(l — f3 3 ). With C — 1, therefore, f3 3 = 1/3 and one obtains the Boldyrev's velocity 
scaling: 

cw „/» + i-ar <8) 



C(3) " V3 

This velocity scaling has been found to provide a very accurate prediction for numerical 
simulations of super-sonic and super- Alfvenic turbulence (Boldyrev et al. 2002a,b). The 
computation of the corresponding structure functions for the gas density distribution is 
discussed in Boldyrev et al. (2002a), where an exponent 77(2) = 0.7 is predicted for the 
second order structure function of the projected density. In this paper we have found the 
'provoking' result that the projected density field of super-sonic turbulence in the ISM has 
a structure function that is almost indistinguishable from Boldyrev's velocity scaling. 

This result should inspire the mathematical work. As shown by Dubrulle (1994) and 
by She & Waymire (1995), a hierarchy of structures producing the scaling relation (5) can 
be obtained by a random multiplicative process with Log-Poisson statistics. The result 
of this paper suggests that the density field of super-sonic turbulence is fragmented by a 
multiplicative process with Log-Poisson statistics. 
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4.1. The (3 and 7 Tests for the Hierarchical Model 

In the previous section the scaling of the structure functions of integrated intensity of 
two maps of the Taurus and Perseus molecular cloud complexes has been compared with 
the velocity scaling of super-sonic turbulence. A different approach would be to use obser- 
vational, experimental or numerical data to determine the parameters (3 and 7 from a fit to 
(5). There is, however, an independent method to first estimate the parameter (3 and then 
obtain 7 from (5), which also tests the validity of two major assumptions of the hierarchical 
structure model at the same time. This method is presented in She et al. (2001); we apply 
it in the following. 

The basic assumption in the derivation of (5) by She & Leveque (1994) is the existence 
of the universal scaling behavior: 

F p+1 (r) = A p F p (rfF (oo) (r) l -P (9) 

where 

F p (r) = S p+1 (r)/S p (r) = (|A/| P+1 )/ (|A/| P ) (10) 

is usually referred to as the pth order intensity of fluctuations and A p are constants inde- 
pendent of r and found to be also independent of p in a number of cases. 

The (3 test verifies the validity of this basic assumption with a log-log plot of F p+ i/F 2 
versus F p /F\. If the plot is a straight line the assumed hierarchical symmetry (9) is satisfied 
and the data pass the f3 test. The plot is shown in Figure 5, top panel, for both Perseus 
(diamond symbols) and Taurus (asterisk symbols). We have combined together values of 
F p (r) for all values of r used to compute the moment scaling exponents r](p) (approximately 
the range of scales between 0.3 and 3 pc). Both Taurus and Perseus data apparently pass 
the (3 test. 

The constants A p are independent of p, as for the hierarchy of intensity of fluctuations 
of velocity in the shell model analyzed by Leveque & She (1997), in the turbulent Couette- 
Taylor flow studied by She et al. (2001) and in laboratory data obtained by Chavarria et al. 
(1995). The slope of the intensity of fluctuation function plotted in Figure 5 provides an 
estimate of the value of (3 = 0.79 for Taurus and (3 = 0.78 for Perseus. This is to be compared 
with the value of f3 — 0.69 of the BSL's velocity scaling for super-sonic turbulence. 

The second basic assumption by She & Leveque (1994) is: 

Foe ~ SI (11) 

which allows them to derive (5). The 7 test verifies directly the validity of (5) by first 
assuming the value of (3 obtained from the (3 test and then by plotting r]{p) — xip, P) versus 
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V — 3x(P)/^), where x(PiP) = (1 ~ /^ p )/(l ~ (3 3 )- If the plot is a straight line the data pass 
the 7 test and the slope of the plot provides an estimate of the value of 7. The plot is shown 
in Figure 5, bottom panel. We are able to fit the plot to a straight line for orders p > 11 
for both Taurus and Perseus. Since 7 is a property of the very high order moments (see 
definition (11)), and since we have obtained a straight line for large moments, we can say 
that the Taurus and Perseus data pass the 7 test. We estimate 7 = 0.11 for Taurus and 
7 = 0.06 for Perseus. 

5. Conclusions 

We have computed the structure functions of the integrated intensity images of two J=l- 
13 CO maps of Taurus and Perseus. The structure functions scale as power laws within the 
range of scales 0.3-3 pc. The scaling exponents have been computed up to the 20th order and 
are statistically significant at least up to the 15th order in Taurus and 8th order in Perseus. 
They are found to follow very closely the velocity scaling of super-sonic turbulence proposed 
by Boldyrev (2002), within 5% and 10% accuracy for Taurus and Perseus respectively. 

We have verified that the projected density field (or integrated intensity) of the Taurus 
and Perseus molecular cloud complexes can be described by a hierarchical model as the one 
proposed by She & Leveque (1994) for the velocity structure functions of incompressible 
turbulence. We have done so by testing the validity of the two basic assumptions of the 
hierarchical model for our data. The validity of the assumptions of the hierarchical model 
means that the integrated intensity images we have analyzed (an approximate estimate of 
the projected density of the Taurus and Perseus molecular cloud complexes) are the result 
of a multiplicative process with Log-Poisson statistics (Dubrulle 1994). 

The complete derivation of the relation between the structure functions of velocity 
and projected density is the subject of future work. However, the close similarity of the 
structure functions of projected density in Taurus and Perseus with that of the velocity field 
of turbulence provides additional evidence that super-sonic turbulence is the major factor 
controlling the density field in the range of densities and scales sampled by the maps we have 
analyzed. 

It is well established that super-sonic turbulence plays an important role in the dynam- 
ics of the cold ISM (Larson 1981; Padoan et al. 1998b, 1999; Padoan et al. 2001; Padoan & 
Nordlund 1999) The statistical properties of this ISM turbulence need to be discovered and 
understood in order to elaborate a statistical theory of star formation (Padoan & Nordlund 
2002) We have shown in the present paper that this can be achieved with existing observa- 



-9 - 



tional data, by studying the projected density field of molecular clouds. This type of study 
will be greatly improved by far infrared imaging of the dust thermal emission from turbulent 
ISM clouds, obtained by future satellite missions such as SIRTF. 

The work of PP was performed while PP held a National Research Council Associateship 
Award at the Jet Propulsion Laboratory, California Institute of Technology. 
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Figure captions: 

Figure 1: Structure functions of orders p = 1 to to p = 20 of the integrated intensity 
J=l-0 13 CO maps of Taurus (Mizuno et al. 1995) and Perseus (Billawala et al. 1997). The 
solid lines are least square fits to the structure functions in the approximate range 0.3-3 pc, 
where the structure functions are well described as power laws. 

Figure 2: Top panel: Structure function scaling exponents relative to the third order, 
Vzip) — r /(p)/ ? 7(3), up to the 20th order, for Taurus (asterisk symbols), Perseus (diamond 
symbols), Boldyrev's velocity scaling for super-sonic turbulence (solid line), She-Leveque's 
velocity scaling for incompressible turbulence (dashed line) and Kolmogorov's velocity scal- 
ing (dotted line). The moments of the structure functions of the integrated intensity images 
of Taurus and Perseus follow very closely Boldyrev's velocity scaling, within 5% accuracy up 
to the 20th order (and beyond) for Taurus and within 10% accuracy for Perseus. Bottom 
panel: Same as top panel, but relative to Boldyrev's scaling, r] 3 (p) BS L- 

Figure 3: Left Panel: Probability density function (PDF) of the integrated intensity 
differences of the Taurus map. Differences are computed at Ar = 2.6 pc. Vertical dotted 
segments show the position of the peak of the integrand of (3) for orders p = 3 to 15. The 
PDF defines well all moments up to the 15th order. Higher moments are converged due to 
the steep cutoff of the PDF. Right panel: Same as left panel, but for the Perseus map and 
for Ar = 2.8 pc. The PDF of Perseus has a very sharp cutoff at large values, which causes 
a fast convergence of high order moments. The PDF of Perseus is such that the peak of 
the integrand of (3) is not very pronounced and a large interval of values of the PDF give 
therefore a significant contribution to the integral. 

Figure 4: Top Panel: The f3 test (see § 4.1) that verifies the validity of the hierarchical 
structure model. Intensity of fluctuations from both Taurus and Perseus are well fit by a 
straight line in this log-log plot. The data therefore pass the (3 test. The value of the pa- 
rameter (3 given by the slope of the least square fits is f3 — 0.79 for Taurus and (3 = 0.78 for 
Perseus. Bottom panel: The 7 test (see § 4.1). The data is well fit by a power law at high 
orders. The solid lines are least square fits for orders within the range p = 12 — 20. The 
data therefore pas the 7 test. The slope of the least square fits provides the values 7 = 0.11 
for Taurus and gamma = 0.06 for Perseus. 
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Fig. 1.- 
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Fig. 2.- 
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Fig. 4.- 



